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Abstract 

This  paper  is  related  to  the  characteristics,  the  development  and  the  simulation  of  a  particular  fluid- 
structure  interaction  phenomenon  -the  continuous  oscillation  of  an  equilateral  triangular  cylinder  in 
the  uniform  incoming  flow.  We  propose  the  explanation  of  this  oscillation.  Preliminary  analysis 
indicates  that  the  cylinder  could  oscillate  incessantly  by  initial  positional  perturbation  or  incoming 
flow  fluctuation.  It  is  the  unbalance  force  acting  on  the  cylinder’s  side  faces  that  causes  such 
movement.  On  one  side,  the  flow  will  be  flow-past-flat-plate  like,  whereas  the  other  side  will  be  flow- 
past-sharp-edge  like.  Due  to  the  unbalanced  pressure  exerting  on  the  two  sides,  the  cylinder  rotates. 
When  the  cylinder  moves,  these  mechanisms  switch  side  interchangeably,  and  bring  the  cylinder  to 
continuous  oscillation. 

Keyword:  fluid-structure  interaction,  oscillation,  triangular  cylinder,  vortex  shedding,  Strouhal 
number,  Finite-element  method 

Introduction 

To  date,  this  coupled  fluid-structure  interaction  has  not  yet  been  well  studied  nor  simulated.  This  is 
due  to  the  unsteadiness  surrounding  flow,  as  well  as,  continuously  moving  boundary  (the  cylinder 
either  oscillates  or  rotates)  around  the  triangle.  The  available  prior  research  works  were  only  dealt  with 
the  stationary  cylinder  at  different  position  (Nakashima,  1992  [ref  1] ),  or  with  other  geometry  (Hu  et 
al,  2002  [ ref  2],  Sakamoto  et  al,  2001  [ref  3],  Naudascher  &  Wang,  1993  [ref  4 ],  Nakamura  & 
Nakashima,  1993  [ref  5]). 

Prior  experiments  (Srigrarom,  1998  [ref  6] )  had  been  conducted  to  investigate  such  behaviour  by 
means  of  Food  coloring  dye  and  Laser  Induced  Fluorescence  flow  visualization.  The  oscillating 
frequencies  appeared  to  be  governed  by  Strouhal  numbers,  which  fell  within  limited  range  of  0.13  < 
Str  <0.18.  Beyond  this  range,  the  cylinder  was  either  stationary  or  rotates  only  in  one  direction 

Due  to  measurement  difficulties  in  experimental  study  of  the  unsteady  flow  pattern  in  the  water  tunnel, 
only  qualitative  results  were  obtained.  Alternatively,  computational  fluid  dynamics  (CFD)  simulation 
presents  a  viable  approach.  Recently,  the  quasi-steady  two-dimensional  flow  study  using  CFD  was 
conducted  recently  to  get  quantitative  results  and  presented  in  this  paper.  The  uniform  incoming  flow 
was  implemented  for  CFD  simulation. 

Proposed  oscillation  model 

Recently,  the  author  has  come  up  with  the  explanation  of  the  self-excited  oscillation  behaviour 
(Srigrarom,  2003  [ref  7])  In  the  following  figure  1,  we  consider  the  flow  pattern  around  a  symmetrical 
cylinder  at  various  positions.  The  upstream  flow  is  uniform  with  zero  angle  of  attack,  and  is  brought  to 
stagnation  at  the  front  of  the  cylinder.  The  flow  is  divided  into  two  (2)  identical  streams,  upper  and 
lower  streams.  The  divided  stream  varies  smoothly  from  a  90-degree  turn  at  S  to  a  stream  merging 
with  the  freestream  at  edge  A  or  B.  At  A  or  5,  the  flow  separates  from  the  cylinder  at  both  tips  or 
shape  edges  of  the  cylinder,  creating  back  flow  or  eddies  on  both  lateral  sides.  At  the  trailing  edge  C, 
the  flow  pattern  changes  according  to  the  freestream  Reynolds’  number.  At  low  Reynolds’  number, 
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Simulation  of  flow  past  an  oscillating  triangular  cylinder  using  Finite  Element  Method 


the  two  divided  streams  create  alternating  vortices,  shed  downstream.  At  higher  Reynolds’  number, 
the  two  streams  join  together  at  alternating  streamwise  locations  creating  a  turbulent  wake 
downstream  of  the  cylinder. 


Fig.  1  (left):  Flow  past  the  equilateral  triangular  cylinder  at  zero  degree  angle  of  attack 
Fig.  2  (right):  Flow  past  the  equilateral  triangular  cylinder  at  positive  angle  of  attack 

In  figures  2  and  3,  supposed  that  the  cylinder’s  symmetry  with  respect  to  the  freestream  is  changed 
slightly.  The  front  face  of  the  cylinder  AB  inclines,  making  an  angle  with  the  freestream.  The  flow  still 
separates  at  both  A  and  B ,  but  the  pattern  is  asymmetrical.  The  radius  of  curvature  of  the  streamline 
separated  at  A,  rA,  is  smaller  than  the  radius  of  curvature  at  B ,  rB.  As  a  result  of  conservation  of  angular 
momentum,  the  velocity  at  A  is  higher  than  at  B\  therefore,  the  pressure  at  A  is  lower  than  at  B 
( pa<pb  ).  The  cylinder  rotates  clockwise  about  the  pivot,  and  the  frontal  surface  AB  becomes  more 
inclined  to  the  freestream. 


Fig.  3:  Flow  past  the  equilateral  triangular  cylinder  at  positive  angle  of  attack,  schematic  diagram 
(continued). 

Due  to  the  continuing  turning,  the  triangle  is  now  is  a  position  where  AC  is  parallel  to  the  freestream. 
The  flow  still  separates  at  the  upper  tip  of  the  cylinder  (5),  but  the  flow  in  the  lower  part,  after  turning 
parallel  to  the  front  face,  separates  only  slightly  at  the  lower  tip  of  the  cylinder  (A),  then  reattaches  to 
the  lower  lateral  face  (AC).  The  lower  flow,  then,  changes  behaviour  to  that  of  a  flat-plate-like  flow. 

By  virtue  of  the  obvious  difference  between  the  two  flow  patterns,  the  local  pressures  differ  at  the 
upper  and  lower  parts  of  the  cylinder.  The  upper  part,  with  the  existence  of  a  large  eddy,  has  lower 
pressure,  especially  at  the  core  of  the  eddy,  compared  with  the  freestream;  whereas  at  the  lower  part 
the  pressure  is  equal  to  the  freestream  pressure,  since  the  fluid  still  flows  smoothly  and  parallel  there. 
Therefore,  the  flow  in  the  lower  part  now  has  greater  pressure  than  the  upper  part,  (PA  >  PB )  and  the 
cylinder  tends  to  rotate  back  to  its  original  position. 
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As  a  consequence  of  the  above  step,  the  unbalanced  pressure  forces  the  cylinder  to  rotate 
counterclockwise  back  to  its  original  position  (under  the  assumption  that  the  cylinder  starts  rotating  in 
a  clockwise  direction,  as  described  in  the  previous  step).  Due  to  the  inertia  of  the  cylinder  and  flow, 
the  motion  of  the  cylinder  does  not  stop  when  it  returns  to  the  original  position  of  figure  1 .  Instead,  the 
cylinder  continues  to  swing  in  the  counterclockwise  rotation.  As  a  result,  the  flow  pattern  becomes  as 
shown  in  figures  4  and  5,  which  correspond  to  those  of  figures  2  and  3  flipped  upward. 


Fig.  4  (left):  Flow  past  the  equilateral  triangular  cylinder  at  negative  angle  of  attack. 

Fig.  5  (right):  Flow  past  the  equilateral  triangular  cylinder  at  negative  angle  of  attack  (continued). 


The  overall  phenomenon  can  be  viewed  as  the  interchange  of  the  flow  patterns,  from  the  flow  past  the 
sharp  edge  to  flow  over  the  flat  plate,  and  vice  versa. 

Computational  Fluid  Dynamics  software:  FEMLAB 

In  general,  the  fluid  dynamics  problems  can  be  described  by  Navier-Stokes  (NS)  equations.  There  are 
a  lot  of  software  programs  available  in  which  NS  equations  are  solved  computationally 
(Computational  Fluid  Dynamics  -CFD).  Most  of  them  are  developed  in  order  to  get  the  simulation 
results  in  a  short  time  with  reasonable  accuracy.  The  disadvantages  are  that  they  are  not  very  flexible 
when  the  user  wants  to  change  the  geometry  or  modify  the  models,  and  they  most  often  act  as  black 
boxes.  Here  we  use  another  category  of  commercially  available  software  like  FEMLAB  [ref.  8]. 
FEMLAB  is  developed  specifically  for  solving  not  only  NS  equations  or  CFD  but  also  other  Partial 
Differential  Equations  (PDEs)  where  the  user  in  principle  can  define  and  simulate  any  system  of 
coupled  PDEs.  Practical  physics  and  engineering  problems  in  the  area  of  fluid  mechanics  and 
structural  mechanics  can  be  solved  with  the  software.  The  problem  solved  in  this  paper  is  quasi-steady 
two-dimensional  flow  of  the  continuous  oscillation  of  an  equilateral  triangular  cylinder  in  the  uniform 
incoming  flow.  One  of  the  main  advantages  of  FEMLAB  is  that  the  user  can  focus  on  the  model  (PDE 
coefficients  on  the  domain  and  the  domain  boundary)  and  does  not  have  to  spend  much  time  on 
solving  and  visualization. 

FEMLAB  is  a  toolbox  written  in  MATLAB  [ ref  9 ].  It  solves  system  of  coupled  PDEs  of  upto  32 
variables.  The  specified  PDEs  may  be  non-linear  (e.g.  NS  equations)  and  time  dependent  and  act  on  a 
one-,  two-  or  three-dimensional  geometry.  The  PDEs  and  boundary  values  can  be  represented  by  two 
forms.  The  coefficient  form  is  as  follows: 


da  —  -  V  •  (cVw  +  ecu  -  /)  +  J3V u  +  au  =  f 
dt 

in  Q 

(la) 

n  ■  (cVu  +  au  -y)+  qu  =  g  -  A 

on  dQ 

(lb) 

hu  =  r 

on  dQ 

(lc) 

The  first  equation  (la)  is  satisfied  inside  domain  Q  and  second  (lb)  -so  called  “generalized  Neumann 
boundary”,  and  third  (lc)  -so  called  “Dirichlet  boundary”  equations  are  both  satisfied  on  the  boundary 
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of  the  domain  <9Q,  n  is  the  outward  unit  normal,  and  is  calculated  internally.  A  is  an  unknown  vector¬ 
valued  function  called  the  “Lagrange  multiplier”.  This  multiplier  is  also  calculated  internally  and  will 
only  be  used  in  the  case  of  mixed  boundary  conditions.  The  coefficients  da,  c,  a,  y,  a,  f  g,  q  and  r 
are  scalars,  vectors,  matrices  or  tensors.  Their  components  can  be  functions  of  the  space,  time  and  the 
solution  u.  For  a  stationary  system  in  coefficient  form  da  =  0  .  Often  c  is  call  the  diffusion  coefficient, 
a  and  /?  are  convection  coefficients,  a  is  the  absorption  coefficient,  and  y  and / are  source  terms. 

The  second  form  of  the  PDEs  and  boundary  conditions  is  the  general  form: 


da-  +  V-T  =  F 

a  dt 

in  Q 

(2a) 

-n-T  =  G  +  A 

on  dQ 

(2b) 

o 

II 

on  dfl 

(2c) 

The  coefficients  T  and  F  can  be  functions  of  the  space,  time,  the  solution  u  and  its  gradient.  The 
components  of  G  and  R  can  be  functions  of  the  space,  time  and  the  solution  u. 

Given  the  geometry  data,  an  initial  finite  element  mesh  is  automatically  generated  by  triangulation  of 
the  domain.  The  mesh  is  used  for  discretization  of  the  PDE  problem  and  can  be  modified  to  improve 
accuracy.  The  geometry,  PDEs  and  boundary  conditions  are  defined  by  a  set  of  fields  similar  to  the 
structure  in  the  language  C.  A  graphical  user  interface  is  used  to  simplify  the  input  of  these  data.  For 
solving  purposes,  FEMLAB  contains  specific  solvers  (like  static,  dynamic,  linear  and  non-linear 
solvers)  for  specific  PDE  problems. 

The  incompressible  fluid  mechanics  in  two-  and  three-dimensions,  is  modeled  within  the 
incompressible  Navier-Stokes  (NS)  application  mode.  The  time-dependent  NS  equations  can  be 
written  as: 

P~jF~  V  *  rjyu  +  (V uf  )+  p(u  •  V)w  +  Vp  =  F  in  Q  (3a) 

V  •  u  =  0  in  Q  (3b) 

In  terms  of  the  FEMLAB  formulation,  the  NS  equations  are  expressed  by  using  the  PDE  coefficients 
according  to  the  following  table.  Notice  that  the  general  form  is  used  in  this  application  mode: 


Coefficient 

Expression 

da 

r 

F 

p 

-rjiyu  +  (Vw)r) 

0 

F  -  p{u  •  V)m  -  V p 
—  V  •  u 

Simulation  procedures  and  results 

The  cylinder  is  initially  positioned  such  that  the  frontal  face  is  at  -30  degrees  angle  to  the  undisturbed 
freestream.  The  body  fit  elements  are  created,  using  the  built-in  function  inside  FEMLAB,  and  refined 
such  that  there  are  always  12,000  elements  in  the  computational  field  to  insure  accuracy.  The  example 
of  the  meshing  (elements)  for  a  particular  wedge  position  is  shown  in  figure  6. 

The  simulation  takes  place  to  get  the  streamline  and  the  pressure  acting  on  the  lateral  faces  of  the 
cylinder  (both  upper  and  lower).  Then,  from  the  pressure  information,  we  can  find  the  moment  acting 
around  the  cylinder  (Srigrarom,  2003  [ ref.  7]).  It  is  the  unbalance  pressures  (moments)  on  the  upper 
and  lower  side  faces  of  the  cylinder  that  lead  to  the  angular  movement  of  the  cylinder  around  the 
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geometric  center.  The  example  pressure  contours  are  shown  in  figures  7  and  8,  for  the  cylinder  at  +30 
degree  and  -30  degree  positions. 


Fig.  6:  Body-fitted  element  around  the  cylinder.  The  total  number  of  element  is  always  around  12,000. 

Pressure  on  the  side  faces  of  the  triangle,  at  +25  deg  position  Pressure  on  the  side  faces  of  the  triangle,  at  -25  deg  position 


Fig.  7  (left)  and  8  (right):  Pressure  acting  on  the  lateral  sides  of  the  cylinder  at  +25  degrees  (left)  and 
-25  degrees  (right)  position.  The  flat  plate  like  and  separation  faces  are  referred  to  their  flow  patterns. 
They  are  either  upper  or  lower  side  faces  according  to  cylinder  position  (see  the  simulation  in  Fig.  9). 


From  the  pressure  unbalance  (and,  hence,  the  turning  moment),  we  can  find  the  consequent  angular 
displacement  at  the  fixed  time  increment.  Then,  we  redraw  the  geometry  -rotate  the  cylinder  to  the 
new  position,  and  do  another  computation  for  streamline  and  pressure  distribution.  The  process  is  then 
repeated  until  the  cylinder  moves  to  another  extreme  position  of  +30  degree,  then,  the  simulation  stops. 
The  termination  of  simulation  at  +30  degree  is  because,  at  this  position,  the  flow  pattern  is  the 
vertically  flipped  image  of  the  initial  position  (-30  degree)  already,  and  that,  one  can  argue  that  the 
reverse  process  will  take  place.  The  enlarged  streamline  flow  patterns  from  -30  degrees  to  +30 
degrees  positions  are  shown  in  figure  9  (a)-(m)  below. 

The  FEMLAB  result  in  figure  9  (a)  shows  that,  on  the  windward  side  of  the  cylinder,  a  vortex  is 
formed  at  the  sharp  edges  of  the  cylinder  during  the  initial  phase  (at  lower  face),  whereas  on  the 
leeward  side,  the  flow  stays  attached  (at  upper  face).  The  size  of  the  separated  vortex  (lower  side) 
increases  with  the  inlet  velocity.  Toward  the  rear  of  the  cylinder  on  the  leeward  side  (upper  side),  the 
sharp  edge  behind  the  cylinder  creates  another  vortex  -of  the  opposite  sign  to  the  windward  one,  and 
causes  reverse  flow.  Behind  the  cylinder,  the  pair  of  windward  and  leeward  vortices  interact  each 
other,  are  shed  downstream  by  the  main  flow,  disintegrated  by  viscosity  and  disappear. 
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(a)  at  -30  degree  (b)  at  -25  degree 


(d)  at  -15  degree  (e)  at  -10  degree  (f)  at  -5  degree 


(g)  at  0  degree  (h)  at  +5  degree  (i)  at  +10  degree 


(j)  at  +15  degree  (k)  at  +20  degree  (1)  at  +25  degree 

Fig.  9  (a)-(m):  The  simulation  of  the  rotating  equilateral  triangular  cylinder  in  clockwise  direction.  The 
Reynolds  number  based  on  cylinder’s  width  =  1. 
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(m)  at  +30  degree 


Fig.  9  (a)-(m)  (Continued):  The  simulation  of  the  rotating  equilateral  triangular  cylinder  in  clockwise 
direction.  The  Reynolds  number  based  on  cylinder’s  width  =  1. 

As  time  goes  on,  as  shown  in  figures  9(e),  (f)  and  (g),  the  flow  on  the  upper  side  becomes  more 
diverging,  and  that  the  new  separating  vortex  starts  to  form.  The  flow  behavior  on  this  upper  side 
changes  from  flat  plate  like  to  separation  like  pattern.  In  the  mean  time,  the  flow  on  the  lower  side 
becomes  more  converging,  and  that  the  existing  separating  vortex  is  suppressed.  The  flow  behavior  on 
this  lower  side  changes  from  separation  like  to  flat  plate  like.  The  flow  on  both  sides  change  their 
patterns,  as  the  cylinder  is  at  zero  degree  position,  as  shown  in  figure  9  (g),  and  even  (h).  Although, 
this  series  show  only  clockwise  motion,  the  interchange  of  flow  behavior  also  exists  at  zero  degree 
position  in  the  returning  counterclockwise  motion,  as  well. 

Insofar  the  simulation  results  have  been  showing  the  self-excited  oscillation  behavior  of  the  cylinder. 

It  is  found  from  the  prior  experiments  that,  the  agreement  extends  not  only  the  overall  oscillating 
characteristics,  but  also,  the  flow  details.  In  figures  10,  the  comparison  of  the  flow  patterns  between 
the  one  from  the  real  flow  phenomena  in  water  and  another  one  from  the  present  CFD  simulation 
(FEMLAB  result),  at  the  same  angular  position  of  +25  degree.  The  shape  of  the  stretched  bubbles  in 
the  water  shown  in  the  left  picture  depicts  the  flow  direction,  and  agree  well  the  streamline  in  the 
simulated  CFD  result  shown  in  the  right  picture.  In  figures  1 1,  the  comparison  of  the  flow  patterns 
between  the  one  from  the  real  flow  phenomena  in  water  and  another  one  from  the  present  CFD 
simulation,  at  the  same  angular  position  of  +25  degree.  The  path  of  the  stretched  bubbles  in  the  water 
shown  in  the  left  picture  depicts  the  flow  separation,  at  the  lower  tip,  and  agree  well  the  streamline  in 
the  simulated  CFD  result  shown  in  the  right  picture.  As  anticipated,  figures  1 1  appear  as  the  vertically- 
flipped  image  of  figures  10. 


fig.  10:  Comparison  of  the  flow  pattern  between  the  real  flow  phenomena  in  water  and  the  present 
CFD  simulation,  at  the  same  angular  position  of  +25  degree.. 
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Fig.  11:  Comparison  of  the  flow  pattern  between  the  real  flow  phenomena  in  water  and  the  present 
CFD  simulation,  at  the  same  angular  position  of  +-25  degree. 


Conclusions 

The  present  CFD  study  of  a  simplified  quasi-steady  flow  past  an  oscillating  triangular  cylinder  found 
an  interesting  self-excited  oscillation  of  the  cylinder  as  observed  in  experimental  investigation  before. 
Under  the  simple  uniform  incoming  flow,  the  interchange  flow  pattern,  at  certain  flow  velocity,  was 
found.  The  results  from  Finite  Element  Method  (FEMLAB  software)  appeared  to  agree  with  the 
proposed  underlying  physical  explanation  of  the  phenomena.  It  is  the  unbalance  force  acting  on  the 
cylinder’s  side  faces  that  causes  such  movement.  On  one  side,  the  flow  will  be  flow-past-flat-plate  like, 
whereas  the  other  side  will  be  flow-past-sharp-edge  like.  When  the  cylinder  rotates,  these  mechanisms 
switch  side  interchangeably,  and  bring  the  cylinder  to  continuous  oscillation.  The  qualitative 
agreement  between  the  theory,  simulation  and  the  experiments  are  pronounced.  This  study  of  the  flow 
past  an  oscillating  triangular  cylinder  could  be  further  enhanced  by  considering  more  robust  fluid-solid 
interaction,  incoming  flow  fluctuation  and  non-Newtonian  effect. 
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